Photoassociation and optical Feshbach resonances in an atomic Bose-Einstein condensate: treatment 

of correlation effects 
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In this paper we formulate the time-dependent many-body theory of photoassociation in an atomic Bose- 
Einstein condensate with realistic interatomic interactions, using and comparing two approximations: the first- 
order cumulant approximation, originally developed by Kohler and Burnett [Phys. Rev. A 65, 033601 (2002)], 
and the reduced pair wave approximation, based on a previous paper [Phys. Rev. A 68 033612 (2003)] general- 
izing to two channels the Cherny-Shanenko approach [Phys. Rev. E 62, 1046 (2000)]. The two approximations 
differ only by the way a pair of condensate atoms is influenced by the mean field at short interatomic separa- 
tions. For these approximations we identify two different regimes of photoassociation: the adiabatic regime and 
the coherent regime. The threshold for the so-called "rogue dissociation" [Phys. Rev. Lett. 88, 090403 (2002)] 
(where mean-field theory breaks down) is found to be different in each regime, which sheds new light on the 
experiment of McKenzie et al [Phys. Rev. Lett. 88, 120403 (2002)] and previous theoretical calculations. 

We then use the two approximations to investigate numerically the effects of rogue dissociation in a sodium 
condensate under conditions similar to the McKenzie et al experiment. We find two different effects: reduc- 
tion of the photoassociation rate at short times, and creation of correlated pairs of atoms, confirming previous 
works. We also observe that the photoassociation line shapes become asymmetric in the first-order cumulant 
approximation, while they remain symmetric in the reduced pair wave approximation, giving the possibility to 
experimentally distinguish between the two approximations. 
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L INTRODUCTION 

The possibility to create molecular condensates opens 
new research avenues. These include test of fundamental 
symmetries QJ 01, determination of fundamental constants 
through molecular spectroscopy with unprecedented accuracy, 
creation of molecular lasers and invention of a coherent 
super-chemistry |3] at ultralow temperatures. Since direct 
laser cooling or sympathetic cooling appear more difficult 
for molecules than for atoms, many experimental groups 
have worked on procedures to transform an atomic quantum 
degenerate gas into a molecular condensate. Over the last 
two years a wealth of new experimental results have appeared 
about formation of molecules in a degenerate gas, starting 
either from an atomic Bose-Einstein condensate fl E l 0. 
or from an atomic Fermi gas j9l fTol. fTTL fT2l fl~3l f 1 4t 1 1 5ll - 

So far the most successful scheme to produce molecules 
in a condensate involves an adiabatic sweep near a Feshbach 
resonance by varying in time the strength of an external 
magnetic field. However, the resulting Feshbach molecules 
are usually in a highly-excited vibrational level of the ground 
state interatomic potential and decay rapidly due to collisional 
quenching. An alternative would consist in varying the 
frequency of a laser in a photoassociation experiment, thereby 
sweeping across an optically induced Feshbach resonance. 
Such resonances have already been discussed in theoretical 
papers (H 113 El and explored in recent experiments 
111 mlH. These resonances are in many ways similar to 
the magnetic Feshbach resonances, except that the resonant 
state is an electronically excited state with a usually very short 
lifetime due to spontaneous emission. As we shall see in this 



paper, this qualitatively changes the many-body properties of 
the system. On the other hand, photoassociation offers more 
experimental possibilities by varying in time the frequency 
or/and the intensity of the laser 1 22] . In particular, the use 
of shaped laser pulses opens the way to better control in this 
domain. For instance, by using a series of shaped pulses it 
should be possible to quickly transfer vibrationally excited 
molecules created via a magnetic Feshbach resonance to their 
ground vibrational state 12311 . Photoassociation with pulsed 
lasers could therefore solve the problem of the short lifetime 
of molecular condensates, in particular bosonic dimers. 

Mixed atomic and molecular condensates, formed by mag- 
netic or optically induced Feshbach resonances, are typically 
described in the mean-field approximation, corresponding to 
two coupled Gross-Pitaevskii equations O. 1241 1251 . This may 
not be sufficient when correlations play a significant role. For 
instance , they must be introduced in the theoretical models 
jp.E7ll2ll29ll30ll to reproduce the damping in the observed 
1 4] oscillations between the atomic and the molecular compo- 
nents of a condensate exposed to a time-dependent magnetic 
field. 

In the case of an isolated resonance, it may be sufficient 
to describe the microscopic quantum dynamics with effective 
interactions such as contact or separable potentials, involv- 
ing parameters fitted on two-body calculations. However, in 
the more general case of photoassociation with shaped laser 
pulses, many levels may be involved; several theoretical studies 
I31LI32II33I1 investigating photoassociation with chirped pulses 
at a two-body level have shown a great sensitivity to the details 
of the molecular potentials. 

A general treatment of photoassociation in a Bose-Einstein 
condensate should therefore be able to take both correlations 
and realistic interactions into account. To achieve this goal, we 
will consider two methods: 
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In a series of papers, |2jJ |3J, |3J, |36j the Oxford 
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group has developed a method based on a truncation 
of the expansion of correlation functions in terms of 
non-commutative cumulants 1 37] . The method to first 
order, hereafter referred to as first-order cumulant ap- 
proximation, has been used so far assuming a separable 
interatomic potential, which was sufficient to success- 
fully interpret magnetic Feshbach experiments 0|35|. 

In a previous paper 1 38], hereafter referred to as paper 
I, we have revisited the treatment of photoassociation 
and Feshbach resonances following another approach, 
hereafter referred to as the effective pair wave approach. 
This work generalizes the ideas of A. Yu Cherny and A. 
A. Shanenko I39l l40ll to two coupled channels. 



The present paper goes further by numerically solving the ef- 
fective pair wave equations as well as the cumulant equations 
with realistic potentials. It is organized as follows: 

• In SectionllTlwe discuss the limitations of contact poten- 
tials and present how to treat realistic interactions using 
either the first-order cumulant or the reduced pair wave 
approximations. 

• In Section |ni| we give the coupled equations for a two- 
channel modelling of photoassociation in a condensate, 
using the two approximations described in SeclITI 

• In Section IIVI we make a connection between these 
equations and the usual mean-field approximation of 
Ref. 12511 for an isolated resonance. This enables us to 
determine the different regimes and conditions for which 
correlations play an important role. 

• Section|V]describes the numerical methods used to solve 
the coupled equations given in Sec. [II]] 

• In Section rVTI we present and discuss the results of our 
numerical calculations for high-intensity photoassocia- 
tion in a sodium condensate, starting from the experi- 
mental conditions of Ref. 1411. 



• We conclude in Section lVTn 

In this paper we give scattering lengths in units of the Bohr 
radius fl = 0.529177 10- 10 m. 



II. DIFFERENT METHODS TO TREAT THE 
CONDENSATE DYNAMICS 

A. Methods with a contact potential 

In many theoretical treatments of the condensed Bose gas 
l42l l43ll . the interaction potential U(r) for two atoms is re- 
placed by the effective contact potential 



4%h 2 c 



■§ 3 (r), 



(1) 



where a is the scattering length of the potential U (r) for parti- 
cles of mass m. The physical argument for this replacement is 



that the detailed structure of the potential U is not resolved at 
the scale of the typical de Broglie wavelength associated with 
very low collision energies. Although the mathematical form 
of the contact potential {1} does not lead to a well defined scat- 
tering problem, it can nevertheless give sensible results within 
theories which treat the interaction perturbatively. For exam- 
ple, in the two-body problem, the Born approximation for the 
scattering length of the contact potential is equal to the scatter- 
ing length a of the original potential 



d 3 rU & {r)=a. 



(2) 



Similarly, in the many-body problem the Hartree-Fock ap- 
proximation for atoms interacting with a contact potential leads 
to the well-known Gross-Pitaevskii equation 1E2I l43ll . How- 
ever, beyond these first-order approximations, ultraviolet di- 
vergences arise due to the zero-range nature of the contact po- 
tential, and one must use a regularized delta function I44l|. or 
renormalize the coupling constant of the contact potential 1 45 ] . 
Once these technical procedures are applied, effective contact 
potentials can be successful in many situations. In particular 
we mention the description of damped oscillations in a mixed 
condensate [27, 28]. 

However, since the effective contact potentials eliminate the 
details of the real potential, some physical quantities, such as 
the kinetic and interaction energies of the gas |40] cannot be 
predicted. Similarly, the presence of bound levels in the po- 
tential is not directly taken into account, whereas these bound 
levels can play an important role in the dynamics of photoas- 
sociation with chirped laser pulses 1 32, 33j]. For these strongly 
time-dependent situations the details of the potentials will be 
important, since their structure is explored over a wide range 
of internuclear separations. We expect that effective contact 
potentials may be inadequate and a realistic interaction poten- 
tial should be used. 



B. Methods using a realistic potential 

Realistic interaction potentials U(r) may have a deep well 
and a strong repulsive wall at short interatomic separations r. 
Because of these features, perturbative treatments like the Born 
approximation are inadequate. For instance, the "Born" scat- 
tering length of the potential U 



a Bom [U} = ^Jd i rU(r) 



(3) 



is a poor approximation of the actual scattering length a. In 
fact, because of the repulsive wall, most potentials are singular 
at r = so they cannot be integrated over space. As a result, 
asorn diverges in principle. In order to obtain finite results we 
need to start from the exact expression for the scattering length 



a[U] 



Anh 2 



d 3 rU{r)<p 2B (r) 



(4) 



where (p?s is the solution of the two-body Schrtidinger equa- 
tion at zero energy: f— + U(r)j q>2s(r) = and is nor- 
malised to be asymptotically equal to 1 . The deviations from 
1 of the wave function ($2B represent the correlations induced 
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by the potential. In particular, because of the strong repulsive 
wall, 92b('") vanishes for r less than a few Bohr radii as it is 
nearly impossible to find two atoms in this region. As a re- 
sult, the integrand U (r)(p2fl(r) in Eq. © is regular for short r, 
which leads to a finite a[U]. Interestingly, Eq. @ can be seen 
as the perturbative expression Q where the "bare" potential is 
replaced by a new effective potential U2B{r) — U (r)(p2e(r). 

In a well-defined many-body theory we also expect that 
when the potential appears in an equation it is always mul- 
tiplied by many-body correlation functions (p^B which go to 
zero for short interatomic separations and are equal to 1 at large 
separations ll39l 14(1 l47il . This introduces effective in-medium 
potentials 1lMB{r) =V '(r)tyMB(r)- At low energy, most proper- 
ties should be determined by the zero-momentum Fourier com- 
ponent of Umb given by the expression 



4%h 2 



d 3 rU(r)y MB {r) 



(5) 



In the region where U(r) is not negligible, many -body ef- 
fects are expected to be small and (p^B should be close to the 
two-body wave function (p?s, so that expression is approxi- 
mately equal to the scattering length a[U] of the two-body the- 
ory. This is consistent with the idea that the two-body scatter- 
ing length determines low-energy properties. 

By analogy with the two-body theory, we may also expect 
that a well-defined many-body theory for realistic interactions 
will have essentially the same structure as a perturbative theory 
with the replacement U (r) — ► Umb^)- This similarity was il- 
lustrated in Refs. ll39l 14811 . where non-perturbative Bogoliubov 
equations were derived and an effective interaction Umb ap- 
peared in place of U. This replacement U(r) — > ZlMB(r) in 
the perturbative theory is equivalent to the contact potential 
replacement U (r) — > t/g(r) as long as we are only concerned 
with low-momentum quantities, such as 0. However, we 
should note that in general these replacements are not equiva- 
lent. Firstly, the Fourier transforms of the effective in-medium 
potentials UMB{ r ) are momentum dependent, whereas for con- 
tact potentials {/g they are constant. Differences might there- 
fore arise for high-momenta, i.e. short distances, the most im- 
portant one being the absence of ultraviolet divergences. Sec- 
ondly, the potentials Umb depend on the many-body dynam- 
ics through ($>mb, whereas the contact potentials t/g are solely 
determined by two-body parameters. In particular, in the time- 
dependent case, the correlation function is time-dependent and 
complex-valued, leading to a time-dependent complex-valued 
effective interaction UMB^t). This feature is not present for 
the contact interaction f/g, which remains constant and real. 
This is another indication that in time-dependent cases, mod- 
els based on a contact interaction or an explicit description 
of the short-range interaction might differ. In any case the 
time-dependent structure of short-range correlations in Bose- 
Einstein condensates is mostly unknown and deserves further 
investigation. 



C. General theory with realistic potentials 

Let us first consider a Bose-condensed atomic gas with 
only one scattering channel. The dynamics of this many-body 
system is given in the second-quantization formalism by the 



Heisenberg equation of motion 

ih^-(x,t) = H x y(x,t) + Jd 3 z^(z,t)U(\z-x\)^i(z,t)^i(x,t), 

(6) 

where \j/(x) is the field operator annihilating an atom at point x, 
H x the one-body Hamiltonian describing the motion of atoms 
in an external potential V(x), 



(7) 



and U is the single-channel interaction potential between two 
atoms (because of the very low density, we neglect interacti- 
ing potentials involving more than two atoms). In the U(l) 
symmetry-breaking picture I43II . the condensate wave function 
*P is given by the quantum average of the field operator 

«P(x,0 = (<\f(x,t)). (8) 

The exact equation for *P is therefore: 

ih—(x,t)=H x x i'(x,t)+ d 3 zU ^|z-x|)(y t (z,f)y(z,f)v(x,0) 

(9) 

According to Eq. ( I80> in the Appendix, assuming that we may 
neglect the influence of the non-condensate atoms colliding 
with the condensate atoms, we can write 

ffl(z,t)y(z,t)y(x,t)) « ^*(z,f)4>(z,x,f), (10) 

where we introduce the quantum average 

4>(x,y,f) = (y(x,?)VKy,?)}- 

This anomalous average is related to the pair wave function 
of two condensate atoms (see Appendix). From Eq. ( 1611 PL we 

get: 

ih—{x,t)=H K *i>{x,t)+ /t/ 3 z1 , *(z,f)f/(|z-x|)4>(z,x,f), 

(ID 

This is a one-body Schrodinger equation with an extra 
term corresponding to the influence of surrounding condensate 
atoms. For convenience, we can express the pair wave function 
<t> in terms of a product of independent particle wave functions, 



*(x,y,f)=Y(vMy,»)<p(*,y,0. 



(12) 



so that all the correlation between two condensate atoms is con- 
centrated into the reduced pair wave function (p l38ll49ll . As- 
suming that <t>(x,y,f) is uncorrected at large distances |x — y|, 
we must have (p(x,y,f) w 1 at such distances. Then we can 
write Eq. (lilt as a Gross-Pitaevskii equation 

Jl", , / 4nh 2 a M (xJ ) , . , . 

/^(x,?) = [H x + ^^|'P(x,f)| 2 )l'(x,f) 

by introducing a time and position dependent scattering length 
aM(x,t), hereafter called the mean-field scattering length, 



a M (x,t) 



4nh 2 



d z C/(|z-x|)(p(z,x,f) 



¥(x,r) 



(13) 
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We can see from this expression that the reduced pair wave 
function cp is the regularising correlation function (|)mb dis- 
cussed in Section 111 Bl In the usual Gross-Pitaevskii theory, 
the pair wave function <t> is assumed to be completely uncor- 
rected: <p(x,y,f) ps 1, which leads to om(x,£) ps asom, where 
a-Born is the Born scattering length defined in Eq. (j3ji. If the 
true scattering length a associated to the potential is close to its 
Born approximation ago™, then we retrieve the standard Gross- 
Pitaevskii equation for V P. However, we have seen above that 
this is not the case for realistic molecular potentials with a re- 
pulsive wall, since asorn(U) — > +°°. We therefore have to keep 
the correlation (p contained in <f>, so that its structure at short 
distance regularizes the singular character of the potential U, 
leading to a value of the mean-field scattering length om close 
to the physical value a. 

The next step is to determine 4>(z,x,f). Starting from the 
Heisenberg equation (jSJ, one can derive the exact equation for 
4>: 

ih— (x,y,t) = (H s +H y + U(\x-y\))<Z>{x,y,t) 

+ S.Jd 3 zU (z-y) (z,f )Vff(z,f )i\r(x,t )ty(y,t )) 

+ {x^y}. (14) 

Eq. (II 41 is a two-body Schrodinger equation with an extra 
term corresponding to the influence of the medium. In or- 
der to take this term into account, we have to determine the 
quantum average (y T (z)y(z)\j/(x)\j/(y)). In fact, Eqs. (tTTl and 
H4\ happen to be the first equations of an infinite set of equa- 
tions coupling each quantum average to higher-order quantum 
averages. Since we want to truncate this set in order to re- 
strict the dynamics to the condensate wave function *P and 
the pair wave function 4>, we need to approximate the average 
(\fr t (z)\j/(z)ij/(x)\}A(y)) in terms of *P and 4>. 

D. First-order cumulant approximation (FOC) 

1. Truncation approximation 

A systematic method to truncate the infinite set of equations, 
named the cumulant method, has been proposed by J. Fricke 
l37ll . This method redefines the different quantum averages 
in terms of quantities called non-commutative cumulants (see 
Appendix), which are constructed in such a way that for a suf- 
ficiently weak interaction, they decrease towards zero whith an 
increasing order. It is therefore natural to neglect, in the infinite 
set of equations, all the cumulants of order larger than a given 
order n. In this way, a consistent truncated set of equations can 
be obtained for any desired order n. 

Doing so, one assumes that the system is close to the ideal 
gas. However, this assumption might fail for realistic interac- 
tions which can induce strong changes in the correlation struc- 
ture at short distances. To extend the validity of the cumulant 
method, T. Kohler and K. Burnett 13411 have devised a modi- 
fied version of the truncation scheme where the free evolution 
of the cumulants of order n + 1 and n + 2 is taken into account. 
This renormalizes the interaction in the first n equations, in a 
way similar to the regularization we discussed in Section lTTBl 



As a consequence, this extended cumulant approach is tech- 
nically applicable to realistic potentials, although it has been 
used so far only with effective separable potentials in the con- 
text of atomic and molecular condensates 1 3 a. 15011. We show 
below how it can be applied to realistic potentials. 

Writing the cumulant equations to the first order, one finds 
Eqs. ( II II and \\A\ with the following approximation: 

ffl(?,t)W*,tW(xW(j,t)) «¥*(v)<s(z,y,*mv)- (15) 

The resulting coupled equations between the condensate 
wave function and the pair wave function can be written as 
follows 160] 



(16) 



m—(x,t) 



(H x +M(x,t))V(x,t) 

(// x +// y + t/(|x-y|))4>(x,y,r)- 
(M(x,t)+M(y,t))V(x,t)¥(y,t) 

involving the mean-field function 

M(x,t) = ^i^x,?)! 2 . 



(17) 



(18) 



These equations have been recently rederived from pair wave 
function considerations |49]. We notice that the mean field M 
acts as an extra potential in Eq. d!6l > but as a source term in 
Eq. Mil . We can verify that in Eqs. ( 1161171 the interaction 
potential U is indeed multiplied by the pair wave function <t>, 
which vanishes at short distances if the correct boundary con- 
ditions are ensured: these equations can therefore be used with 
a realistic molecular potential. 



2. Stationary states 

The cumulant method is designed mainly to propagate the 
dynamical equations during a limited period of time I l34ll . It 
is not well designed for stationary states, which would require 
that the equations are valid for an infinite period of time. How- 
ever, in order to solve the time-dependent equations ( 1161171 . 
an initial condition is required. Here, one cannot assume that 
the pair wave function is initially uncorrelated, as is proposed 
in Ref. 13611 : in the case of realistic interactions such a choice 
would lead to divergence, as we explained in Section UlCl If 
we start from equilibrium, the initial condition has to be cho- 
sen as close as possible to a physical stationary state, typically 
given by the stationary Gross-Pitaevskii equation. 

Setting *P(x,?) = x t , (x)e _!A " ///i , where /j is identified with the 
chemical potential we have to set <t>(x,y.f) = '&(x,y)e~ l2fM l h 
accordingly to obtain, from Eqs. ( I16I17I I, stationary equations 
for and <i>. Introducing M(x) = M(x,0), these equations 
read: 

^(x) = (// X +M(x))*(x) (19) 
2^(x,y) = (// x + // y + C/(|x-y|))<i>(x,y) 

+ (M(x)+M(y))*(x)*(y), (20) 
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and, as explained in Ref. l49Tl . lead to the stationary Gross- 
Pitaevskii equation for /j — ► 0. However, in typical experi- 
ments, fj can be much larger than the energy level spacings 
of the trapping potential V(x). In the experiment of Ref. |41] 
for instance, the trap frequency is about 200 Hz, while the 
chemical potential /j w ° p c (where p c is the density at the 
center of the trap) is about 6 kHz. For chemical potentials 
of this order of magnitude, there is a problem with the solu- 
tion of Eqs. (I19I20I . Indeed, when Eq. ( 1191 is solved for- 
mally using the Green's function G jU (x',y',x,y) of the operator 
H s +H y + U(\x-y\)-2fi: 

4>(x,y) = J d 3 x'd 3 y'G,(x'y,x,y) 

(M(x')+M(y')^(x')V(y'), (21) 

a singularity occurs in this Green's function each time /j is 
close to an energy level of the trap. When Eqs. ( I19I20> and 
i2l\ are solved consistently, this leads to series of unphysical 
resonances of the mean-field scattering length as fj sweeps the 
series of energy levels of the trap (or equivalently, as the den- 
sity of the system is varied). As noted in Ref. |49], this prob- 
lem can be fixed in the limit of a homogeneous system: the 
Green's function becomes undefined for any value of /j, since 
singularities form a continuous spectrum in this limit, and one 
can regularize it either by taking the principal value of the in- 
tegral in Eq. i2li defining the formal solution, or by evaluat- 
ing it slightly outside the singularity line at fj ± i£. The dif- 
ficulty remains however in a numerical implementation, since 
the calculations have to be implemented within a limited box 
(see Section [V), leading to a discretization of the continuum. 
In between the spurious resonances, the mean-field scattering 
length is equal to the correct scattering length a, as we shall 
show in SectionlVl 



Moreover, one can deduce l36ll the very interesting relation 

R'(x,y,t) = Jd 3 z<i>'*(x,z,t)<t>'(z,y,t), (24) 

which, in principle, guarantees the positivity of the non- 
condensate density matrix, so that all occupation numbers of 
the non-condensate modes remain positive. As noted in 1 49], 
this relation is an approximation of the Bogoliubov relation 

fl'(x,y,f)+ Jd 3 zR'(x,z,t)R'(z,y,t) = 

jd 3 z&*(x,z,t)<l>'(z,y,t) (25) 

for small values of the condensate depletion. 

However, we should remark that the initial non-condensate 
fraction predicted by Eq. (1241 is generally an overestimation, 
because of the inappropriate long-range behaviour of <t> within 
the first-order approximation. For instance, in a homogeneous 
system, <t> 1 — a/r at large distances r (which is known to be 
incorrect from the Bogoliubov theory - see Ref. H), so that 
the right-hand side of Eq. (I24> diverges. When calculations are 
performed within a box, the results depend on the size of the 
box in an unphysical way. Therefore one should discard the ini- 
tial number of non-condensate atoms as meaningless, consid- 
ering only the subsequent variations. Note however that these 
variations may become negative, which would be interpreted 
as negative occupation numbers. 



E. Reduced pair wave approximation (RPW) 

1. Approximation on quantum averages 



3. Conservation laws 



Although we have neglected the influence of the non- 
condensate atoms in Eqs. ( I16I17> . these atoms are present in 
the system and described by the second-order cumulant 



R'(x,y,t) = (^(x,t)q(y,t))-V*(x,t)V(y,t), 

which we shall call the non-condensate density matrix. To first 
order in the cumulant truncation scheme, this density matrix 
obeys the following equation: 

dR f 

ih—{x,y,t) = [H y -H x }R , (x,y,t) + 

Jd 3 z <&'*(x,v)t/(|z-y|)<J>(z,y,*) - {x «-> y}*, (22) 

where 4>'(x,y,f) = 0(x,y,f) — *¥(x,t ) l P(y,f) is the correlated 
part of the pair wave function. We can check 1 36] that Eqs. d!6i 
and 1221 lead to the conservation of the total number N of par- 
ticles, 

^J^M^ + J^iW). (23) 



We now start from another approach to close the dynamics 
of Eqs. il 1> and J14i . The pair wave function approach, devel- 
oped in paper I [38], and recalled in the Appendix, shows that 
we can express some quantum averages in terms of eigenvec- 
tors of the two-body density matrix, called pair wave functions. 
If we retain only the contribution from the condensate atoms, 
(thus neglecting that of the non-condensate atoms), we can ex- 
press the quantum averages of interest in terms of *P and <t> 
exclusively. For instance, the quantum average which appears 
in Eq. fl!4i can be approximated by 



0(z,x,f) < l>(z,y,f)*(x,y,f) 



>F(z,f) 



V(x,t)V(y,t)y(z,t) 



(26) 



Compared to Eq. dl5> . this approximation contains extra cor- 
relations between the coordinates, in a way which respects the 
symmetry of the quantum average by permutation of the last 
three coordinates. Note that in Eq. (I14i this average is multi- 
plied by the interaction potential U(\z — y|). Moreover, since 
Eq. (I14i is the equation of motion for the pair wave function 
<I>(x, y), we suspect that the correlation between the positions z 
and y, or between the positions x and y, are the most important 
ones. Neglecting the correlation contained in <t>(z,x), we get 
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the expression 

(¥ + (v)¥(v)¥(x,f)¥(y,0) ~ 



*(z,y,0*(x,y,0> (27) 



which now breaks the symmetry. 

The resulting equations can then be written 

8>P 



(x,r) = (i/ x +M(x,f)) v I , (x,f), (28) 

34>, 



/£— (x,y,f) = [H x +H y + U(x-y) + 

M(x,r)+M(y,f))4>(x,y,r), (29) 

where the mean field M(x,t ), defined as in Eq. (II 81 . now plays 
the role of a potential in both equations. Comparing with 
Eq. i!7i . one sees that a pair of condensate atoms "feels" the 
mean field as a potential not only at large interatomic distances, 
but also at shorter distances where the two atoms are corre- 
lated. This is the only difference with the first-order cumulant 
equations ( 1161171 : although it may seem minimal, we shall see 
in Section lVTl that it may have noticeable consequences on the 
dynamics. 

Using Eq. (1 1 4i and the approximation 127/1 . we find the equa- 
tion for (p: 



r, 2 



2m 



(V2 + V2 



Cy + C x - 



t/(|x-y|)j(p(x,y,f), 
(30) 

with C y — ^^f V y In *P(y ) ■ Vj, . Since the condensate wave func- 
tion VP is expected to be almost uniform in the domain of vari- 
ation of cp, we can neglect the crossed terms C y and C x in 
Eq. (I30i . This condition is of course exact in a uniform system, 
where the crossed terms are stricly zero. In such conditions, 
there is no momentum transfer between the wave function of 
the condensate and the reduced pair wave function: Eq. (1301 is 
then simply the two-body Schrodinger equation in free space. 
In principle, this equation should contain other terms describ- 
ing the influence of the medium 1 38 ] : since they are neglected 
within the approximation ( 1271 . we call it the reduced pair wave 
approximation. 



2. Stationary states 

The stationary states are easily found. Setting again 
*F(x,0 = ¥(x)e-W H and 4>(x,y,f) = 4>(x,y)e- , ' 2 '"/ ri , which 
implies <p(x,y,f) = <p(x,y,0) = cp(x,y), we find: 

f&(x) = (# x +M(x))<P(x), (31) 
= [-^(v2 + v2) + f/(|x- y |)]cp(x, y ) (32) 

From 1321 . we see that independently of the value of 
the chemical potential p, <p(x,y) is the stationary solution 
92b(|x — y|) of the two-body scattering problem at zero energy, 
with the limiting condition (p — > 1 at large distances. As a re- 



sult, the mean-field scattering length 1131 1 is given by: 



5m(x) = ^y^ (|z - x|)92s(|z - x|) H# (33) 



4kH 2 



d 3 rU(r)<p 2B (r) 



(34) 



where we assumed, as discussed in paper I |38], that the con- 
densate wave function does not vary significantly at the scale 
of the potential range. According to Eq. (0}, the right hand 
side of Eq. J34I is the scattering length a, so that Eq. 131> 
is nothing but the usual stationary Gross-Pitaevskii equation. 
To this respect, the reduced pair wave approximation 127> is 
free of the difficulty encountered with the first-order cumu- 
lant approximation 115t . This improvement is obtained be- 
cause the truncation of the higher-order quantum averages is 
less strict: the reduced pair wave approximation includes the 
extra term (M(x) +M(y))<I>'(x, y ) coming from the three-body 
wave function. Note that this term is also present in its pertur- 
bative form jd 3 Z (U(\z-x\) + U(\z-y\)) |*P(z)| 2 <I>'(x, y ) in the 
Hartree-Fock-Bogoliubov approximation. This suggests that 
a more consistent extension of the reduced pair wave approx- 
imation should generalize all the terms of the Hartree-Fock- 
Bogoliubov equations to nonpertubative expressions. 



3. Conservation laws 

Using Eqs. (|6}, we can deduce the exact equation of motion 
for the non-condensate density matrix: 

dR' 

ih—(x,y,t) = [H y -H x }R'{x, yi t) + 

d\ U(\z-y\)(^ f {x,t)^{z,t)^f{z,t)^{j,t)) 

-vp*(x,0( V t( z ,i) ¥ ( z ,i) ¥ (y,/))) | - {x ^ y}*. 

Retaining only the contribution from the condensate atoms, 
Eqs. (18018 li from the Appendix show that R' satisfies Eq. 122> . 
so that the conservation equation (I23i still holds. In con- 
trast, the approximate Bogoliubov relation 124t is not fulfilled 
any longer. We have not been able so far to find an alterna- 
tive relation which would guarantee the positivity of the non- 
condensate occupation numbers. 



III. TWO-CHANNEL COUPLED EQUATIONS FOR 
PHOTOASSOCIATION 



Coupled equations for two atoms with a photoassociation 
laser 



Let us now turn to the standard two-channel model describ- 
ing the photoassociation reaction in a thermal gas of atoms, 
schematized in Fig.^ with a ground potential U g (r), an excited 
potential U e (r) +E e and a time-dependent radiative coupling 
due to a laser red-detuned by HA relative to an atomic reso- 
nance line. We assume that t/ e (°°) — > E e , C/g(°°) — > 0, E e being 
the dissociation limit of the excited potential, while the origin 
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Ground-state 
channel (g) 



Figure 1 : Scheme of the photoassociation process 



of energy is the dissociation limit of the ground potential. The 
relative motion of two atoms is described by a two-component 
wave function <t>, with F gmun i(r,t) describing two atoms collid- 
ing in the potential U s (open channel), and F exc (r,t ) describing 
the ro-vibrational motion in the excited potential U e (closed 
channel). For a continuous laser with constant frequency co 
such that h(0 = h(0o — A, in a classical model, the electric field 
oscillates as cos(cof) and is coupled with the transition dipole 
moment D(r) between the two electronic states. At low tem- 
perature, only i-wave scattering should be considered, and we 
shall neglect rotation effects and transfer of angular momentum 
between the light field and the pair of atoms, so that instead of 
(^ g round(r,f), ^exc(r,f)) only radial functions <j) g (r,f), $ exc (r,t) 
need to be introduced. Using the rotating wave approximation 
15 ill , it is possible to eliminate the rapid oscillations in the cou- 
pling term by defining new wave functions 



$ e (r,t) =exp(/cof)(j)exc('-,0 



(35) 



so that the the time-dependent Schrodinger equation becomes 

in 



dt V <t>e( r >0 



(36) 



t + U g (r) W(r,t) \ (h(r,t) 

W(r,t) t + U e { r )+E e -h®)\$e{r,t) 



where we have introduced the kinetic energy operator T and 
a coupling term: 



B. Coupled equations for a condensate 

Considering now photoassociation in a condensate, we treat 
the many -body problem with a two-component pair wave func- 
tion <t>: the usual component <f> g in the ground (open) chan- 
nel, and a new component <J> e in the excited (closed) channel. 
The component in the open channel is connected at large dis- 
tances to a product of condensate wave functions whereas 
the component in the closed channel corresponds to a bound 
vibrational level and vanishes at large distances, 



* g (x,y,0 
*e(x,y,f) 



l x -y|- 







For simplicity, we will consider a homogeneous system. In 
this case, the condensate wave function *P(x,?) is uniform, and 
the functions <I>(x,y,f) and R'(x,y,t) depend only on the rela- 
tive coordinate x — y = r. Assuming an isotropic situation with 
s-wave scattering only we shall simply write *f?(t), <P g (r,t), 
<t> e (r,t) and R'(r,t) for the condensate wave function, the two 
components of the pair wave function, and the density matrix 
for the non-condensate atoms in the ground state. 



1. Coupled equations in the first-order cumulant approximation 



The first-order cumulant equations (I16I17> are generalized 
as follows [50]: 



dt 

lR dt U e (r,0 



with the mean field: 



- 2M(f ) 



(39) 



(40) 



M(t) = ^Jdh¥*(t)(u e (r)<P g (r,t)+W(r,t)Mr,t) 

(41) 

now depending upon the radiative coupling besides the ground 
state potential, and upon both components of the pair wave 
function. From M(t) one can deduce a mean-field scattering 
length through 



W(r,t)nW(t) = ~J^D. 

2 V ce 



(37) 



I(t) is the intensity of the laser, the constants c and £o are re- 
spectively the ligth velocity and the vacuum permittivity, and 
D is one component of the dipole moment operator, depend- 
ing upon the polarisation of the laser, and assumed to be r- 
indedependent. Since E e — fi,co=A, the 2-channel Hamiltonian 
in (I36> reads 



{r,T) \W(t) t+U e (r)+Aj 



(38) 



The relevant quantities for the photoassociation reaction are 
illustrated in Fig.Q 



a M {t) =M(f)- 



' Anh 2 \V>{t)\ 2 ' 
The conservation equation reads 

p = P g (f) + pg(f) + Pe(f) 



(42) 



(43) 



where p is the total density, p a (f) = ^(f)! 2 is tne density of 
condensate atoms in the initial channel, p' s (t) = R'(0,t) is the 
density of non-condensate atoms in the initial channel, and 
p e = 4-n [\<t>e(r,t)\ 2 r 2 dr is the density of atoms in the excited 
(molecular) channel. Furthermore, we still have in the ground 
channel the "approximate Bogoliubov relation" 1241 . which 
now reads: 



R'(r,t)= ld 3 z&*(z,t)& g (r-z,t). 



(44) 
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and links the density matrix for the non-condensate atoms to 
the correlated part of the component of the pair wave function 
in the ground channel. 



2. Coupled equations in the reduced pair wave approximation 
The generalized form of Eqs. ( 1281291 is 

ih^l = M(tW(t) (45) 
at 



(46) 



with the mean field M(t) still defined by Eq. (I41> . The con- 
servation equation ( 1431 still holds, but the Bogoliubov relation 
J44i is no longer valid. 



V / 2 l P m (0 ( P;)j('")' where cp m (r) is the volume-normalized wave 
function for the relative motion in this resonant molecular level 
and *¥ m (t) is the (here uniform) centre-of-mass wave function 
corresponding to a molecular condensate. With this normaliza- 
tion, | 'I',,, | 2 is the density of molecules. 

Then, the crucial point is to write the correlation in the pair 
wave function as a sum of an adiabatic correlation and a dy- 
namic correlation 

<t> g (r,t) = *P 2 (r)+4>f (r,*) +**"(/•,*) 

= ^ 2 (f)(l+(pf(r,f)+(pf"('-,0)- 



The adiabatic part is found by setting <t>g y " and g^™" to zero 

in Eq. d40l or tyf" 1 and ^tff to zero in Eq. d46l . We can then 
expand the correlation in terms of the scattering states |<pk) of 
the ground potential U g : 



C. Introduction of the spontaneous emission 

When a pair of atoms in a confining trap is photoassociated, 
populating a bound vibrational level of the potential U e (r), the 
excited molecule has a finite lifetime and decays back to the 
ground electronic state by spontaneous emission of a photon. 
Most often, the final state is a continuum level of the ground 
potential U g (r), where the pair of atoms has enough energy 
to escape the trap. In some cases, the radiative transition 
populates a bound level of U g (r), leading to the formation of 
stable molecules I52,l53ll. 



In the previous equations, decay by spontaneous emission is 
not considered. Assuming that it is mainly a loss phenomenon, 
it can be accounted for by adding an imaginary term to the 
excited potential U e (r), 



U e {r) 



(47) 



where 2k /j is the radiative lifetime of the bound levels in the 
excited potential (we have assumed the dipole moment to be r- 
independent). The component of the pair wave function in the 
closed channel therefore contains an exponentially decreasing 
factor: as a result, the total density J43i is no longer conserved 
during the time evolution. 



IV. MEAN-FIELD EQUATIONS AND ROGUE 
DISSOCIATION 



with the normalization (cp k | cp q ) = (27t) 3 8 3 (k-q) l6lll . We find: 



(48) 



Eliminating this adiabatic part in Eqs. J39I40I or d45l40i leads 
to the terms of the mean-field equations 1 38 ] . The equations 
for *P and then read: 



go^ 2 



d 3 k 



im m = hid'-i^yn+^Mttt* 



75iF gkC k 
2 



<tl\n 



(49) 



where 



gk = J d 3 rU (r)cpk(r) is the atom-atom scattering cou- 



4nh- 



pling constant (in particular go 

• Wk = \/2(<Pk|W|<Pm) is the atom-molecule coupling con- 
stant 

• 8' = 8 + E se if is the detuning shifted by the self-energy 
of the molecules 

The equation for the dynamic correlation in the first-order cu- 
mulant approximation is: 



i fit,: — 



fi 2 k 2 A 



(51) 



As a general rule, the coupled mean-field (Gross-Pitaevski!) 
equations of Ref. 12511 are retrieved from both approximations 
when the pair dynamics can be eliminated adiabatically with 
respect to the one-body dynamics 1 38 ] . However the adiabatic- 
ity condition given in Paper 1 13 811 (merely comparing the cou- 
pling constants in the equations) has to be refined. We give 
here a more relevant condition, related to the rogue dissocia- 
tion analysis of ^Ref. 1 54] . 

Let us first assume that only one excited molecular level is 
resonant. In this case, we call 8 the detuning of the laser from 
this molecular level (see Fig. [Q. One can write <t> e (r,f) w 



and in the reduced pair wave approximation: 



. t .dyn fa k dyn -t- -ad 

i he.' = c. —incvr 

k m k k 



(52) 



with ct = Ck/ 1 ? 2 . As long as the momentum distribution of 
C k ' 1M lies in the Wigner's threshold regime, we can make the 
simplification: gk sa gQ and Wk ~ wo. Note that this approx- 
imation does not lead to any ultraviolet divergence, because 
we have eliminated the adiabatic correlation beforehand. Such 
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an approximation, reminiscent of the use of a contact potential, 
would have given rise to ultraviolet divergences if we had made 

lor 11451401 . 

that when the dynamic 



it in the original equations ( 1391 
It is now clear from Eqs. J49I50I > 



correlation C^r n is negligible, one gets the Gross-Pitaevskii 



equations of Ref . 12511 : 



(53) 
(54) 



im> m = ft(8'-^)¥„ 



The mean-field approximation will break down when the dy- 
namic correlation is not negligible any more. This means that 
molecules break into pairs of atoms which will contribute to 
the non-condensate fraction instead of the condensate fraction. 
From Eqs. 



J49I50I > we see that such rogue dissociation |54] 
can be neglected as long as: 



d 3 k f^dyn 
(27?^ U k 



(55) 



To give a more explicit condition, we first have to make a dis- 
tinction between two different regimes of Eqs. J53l54b . 



A. The adiabatic regime 



We consider the limit when |/fi. v P, n | <C |jWo l I ,2 |. In this case, 
•Prn can be eliminated adiabatically, so that: 



\51\ , we find that the condition l/W,,,! <C | ^ wo x i /2 1 is satisfied 
whenever: 



(60) 



^|8' + ^|»woyp, 2gp 



This corresponds to the off-resonant regime of Ref. 12511 : the 
detuning and the spontaneous emission are large with respect 
to the coherent couplings. 



B. The coherent regime 

On the other hand, we may consider the limit \h(5' — 
< |iw l P 2 | a nd lgo|q f ^1 < jwoW*^]. In this case, 
the solutions of Eqs. J53I54I for an initially all atomic system 
are: 

^ w % /p/cosh(w 0> /p0 
¥ m « -/Vptanh(w 0v / P0 

There is a coherent conversion between atoms and 
molecules. In this regime, *P m scales as ^/p and plays the role 
of a one-body field, i.e. a molecular condensate in its own right. 
As a result, the system cannot be described by a two-body the- 
ory any more. For instance, the pair wave function is still of 
the form (1591 . but the mean-field scattering length % is now a 
time-dependent quantity a — /Lsinh(2wo v / pf) where: 



L = 



m wo 



and: 



U'() 



2h(8'-il) 



(56) 



2h{h'-i\) 



| V |2 v =4jt^oa£ 2 ^ (5?) 
m 



In this regime, the coupled equations reduce to a single Gross- 
Pitaevskii equation |25, 28] where the mean-field scattering 
length au corresponds to the modified scattering length: 



A = a- 



^Opt " 



ib 



Opt 



given by the two-body theory of the resonance (30. There- 
fore, most physical properties can be described by the usual 
two-body theory. For instance, Eq. d57t implies that the atomic 
density p ? = I^Fj 2 follows a simple rate equation: 



opt p 2 



(58) 



and using J48i and J56b . one can calculate the pair wave func- 
tion and find the usual asymptotic behaviour: 



is a many-body length which depends on the density. This 
means that the density of atoms does not follow a rate equa- 
tion. 

We find that this coherent regime occurs when: 



wo % /p>fi|5' + ^| 



2gp 



(61) 



In practice, while it is possible to set the detuning to zero by 
varying the laser frequency, the spontaneous decay y is fixed 
by the molecular level of the system. Spontaneous emission in 
alkali systems, and loss processes in general, make it very dif- 
ficult to reach the coherence condition 16 H with realistic den- 
sities and laser intensites. For this reaso n, p hotoassociation 
experiments in condensates so far Il20l.l4ll. 15511 have been con- 
fined to the adiabatic regime and can be described simply 
by two-body theories. 

However, we made some estimates which indicate that pho- 
toassociation for weakly-allowed transitions, such as those 
found in alkaline-earth dimers, would lead to the coherent 
regime. The reason is that spontaneous emission scales as the 
square of the dipole moment whereas the coherent coupling 
wo^/p only scales as the dipole moment. 



•>P 2 (0(l 



(59) 



C. Rogue dissociation 



As a matter of fact, the molecular field scales as p in 
this regime: it is merely a two-body field playing the role of an 
intermediate state during the collision process. Using ( I56l l and 



In the coherent regime, which is the regime originally inves- 
tigated in Ref. 1 54], the relevant energy scale set by the dynam- 
ics is h£l = WQy/p defining the coherent Rabi frequency Q./2K. 
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Figure 2: Different regimes of photoassociation in a condensate, as a function of laser intensity and the combined effect of detuning and 
spontaneous emission. These graphs correspond to the sodium resonance described in I VI Al for a density of 4 x 10 14 at/cm 3 , but spontaneous 
emission is taken as a free parameter. Left panel: the solid line is the limit between the adiabatic regime and the coherent regime of the mean-field 
equations J53I54I . The mean-field approximation is valid in the white areas, and rogue dissociation occurs in the shaded areas. The vertical 
dashed lines corresponds to the left hand side of 1621 being equal to n 2 and the oblique dashed line corresponds to the left hand side of J63I 
being equal to Tt 2 . The arrows indicate how these lines move when the density is increased. Right panel: larger view showing the experimental 
investigation of Ref. l4lll (from dot a to dot b) where the spontaneous emission parameter y/2 is equal to 2jix9.18 MHz; dot c corresponds to 
the condition of our numerical investigation described in lVI CI 



As a result, Eq. ( 15 U or J52I can be made dimensionless us- 
ing the characteristic time X = Q.~ l and length cj = J h/mQ. = 



(4nLp)-^ 2 . Since ^ - p and leg- 
ation condition (I55> can be written: 



dyn i 



1, the rogue dissoci- 



ha 



3/2 



= ((47tL) 3 p) 1/2 «27t : 



(62) 



This condition was presented in J. Javanainen and M. 
Mackie's work |54]. 

Similarly, in the adiabatic regime, the relevant energy scale 
is for = WqP/^/8' 2 + (y/2) 2 corresponding to the mean field 
energy hT = 4nh J A \ p, Eq. fsTi or d52l can be made di- 
mensio nless u sing the characteristic time To = T 1 and length 
cjo = y/h/mT = (4n\A\p)~ 1 / 2 (note that this is the usual heal- 
ing length), and the rogue dissociation condition is now: 



HT 



3/2 



((47 C jA|) 3 p 



1/2 



<27t 2 



(63) 



Both conditions J62b and d63l > can be written as J pa 3 ^ -C 1 



which is a straightforward generalization of the usual condition 
\J pa 3 <C 1 for the validity of the Gross-Pitaevskit equation. 
This can be interpreted as follows: rogue dissociation occurs 
when the average spacing p~'/ 3 between the atoms becomes 
of the order of the typical length associated to the adiabatic 
correlation, either A in the adiabatic regime, or L in the coher- 
ent regime. However, there is a fondamental difference. A is 
the modified scattering length, a two-body quantity indepen- 
dent of the density; in the adiabatic regime, one can thus reach 



the rogue dissociation regime by increasing the density, so that 
the average spacing between the atoms becomes of the order 
of A. The mean-field equations are therefore valid at low den- 
sity in the adiabatic regime. On the other hand, in the coherent 
regime, the many-body length L decreases with density, and 
more rapidly than the average spacing between the atoms. As 
a result, one should decrease the density to observe rogue dis- 
sociation. The mean-field equations are therefore valid at high 
density in the coherent regime. An overview of the different 
regimes and conditions is given in Fig. |2] 



V. NUMERICAL METHODS AND CALCULATIONS 

We now turn to the numerical resolution of the time- 
dependent equations J39I40I and (1451461 . The main difficulty 
is that they describe simultaneously the microscopic dynamics, 
with short characteristic times, and the macroscopic dynamics, 
with longer ones. At each time step, the mean field M(t) and 
the mean-field scattering length am{t) must be determined 
through Eqs. (I41I42> . from the two components of the pair 
wave function <J> g (r,f) and 4> e (r,f) and from the condensate 
wave function l'(f). Then the mean field influences the 
evolution of the three wave functions in the coupled equations 
(1391401 or (I45I46I . Note that in previous calculations based on 
the cumulant equations 1 50], the use of a non-local, separable 
potential made it possible to first perform an independent 
integration of the two-body equation J46l >. then inject the 
results into the condensate equation to solve a non-Markovian 
non-linear Schrtidinger equation. This simplifying procedure 
is not implemented in the present paper where we use the same 



realistic potentials in both approaches. 

The pair wave functions are represented on a grid, using a 
mapping procedure where the grid steps are adjusted to the lo- 
cal de Broglie wavelengths 15611 in the two potentials, choosing 
the smallest of the two values at a given r. The short range 
oscillations are thus described with a dense grid and the long 
range behaviour with a diffuse one. Typical value of the grid 
length L is 200,000 oq. The whole set of equations, either 
J39I40I > or J45I46I are then solved, using the standard Crank- 
Nicholson method 15711 for the propagation in time. The solu- 
tions at t — are chosen as stationary solutions for the open 
(ground) channel, when no laser coupling is present. 

A. Initial state in the first-order cumulant approximation 

Stationary solutions (in the rotating frame defined by 
Eq. ( I35» of Eqs. J39I40I are given by 

* = y/fi (64) 

(it'O - -(■*<*>-*.)"'■*.(&) <*> 

where p g is the initial condensate density, and /j = M is the 
chemical potential (equal to the initial mean field in the present 
case). When the laser is initially off, the coupling term W 
in H' 2 ) is zero. To compute 1641651 . we start from a given 
value for the condensate density p g and a guess value for /j. 
We determine the pair wave function from d65l > by performing 
a matrix inversion. This gives a new value of the mean field 
through Eq. I41i . and the procedure is iterated until conver- 
gence is reached for /j. As explained in Section [I]] for certain 
values of the density p g , the matrix is singular and the solu- 
tion shows an unphysical resonance of the mean-field scatter- 
ing length. Between the resonances the mean-field scattering 
length is essentially equal to the expected scattering length of 
the two-atom system. This difficulty is illustrated in Fig. [3] We 
therefore choose an initial density that is both close to the ex- 
perimental value and such that the mean-field scattering length 
does lie in-between the resonances. 



B. Initial state in the reduced pair wave approximation 

Stationary solutions of Eqs. 1451461 are obtained by solving 
* = y/Vg (66) 

° - (H (2 '(,.)-(^)).(!gj) m 

The two-component function (<J>g,<i>e) is therefore the zero 
energy stationary solution of the two-channel problem, the po- 
tential U e being shifted down by twice the chemical potential. 
This means that the energy of the optically induced Feshbach 
resonance is shifted by the mean-field energy of the conden- 
sate atoms, even though the collision energy is nearly zero. To 
determine eigenfunctions of Eq. J67> we propagate Eq. ( 146 > in 
imaginary time with the Crank-Nicholson scheme 15711 . The fi- 
nite value of the time step value acts as an energy filter, which 
selects the zero energy scattering state in a few steps. 
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Figure 3: Mean-field scattering length as defined from Eqs. <4 11421 . 
calculated from the stationary solution of the first-order cumulant 
equations < 1611 171 . as a function of the mean-field energy. The numer- 
ical calculations shown here correspond to the case of a sodium con- 
densate without any photoassociation laser, and are performed within 
a square box of width d= 10,000 ao (solid lines) or <i=5,000 ao (dashed 
line). The position of the first two levels of the wider box are indicated 
by arrows on the upper horizontal axis, the second one coinciding 
with the first level of the d=5,000 box. As explained in the text, the 
mean-field scattering length has an unphysical divergence each time 
the mean-field energy coincides with an energy level of the box. The 
latter depend on the arbitrary size the box; in our calculations, this 
size can typically go up to 200,000 atomic units, which increases the 
number of resonances. The typical mean-field energy in the experi- 
ment of Ref. 1 41] is about 6 kHz ~ 300 nK, and indicated by an arrow 
on the lower horizontal axis. 



C. Accuracy check: computation of the scattering length 

Checks on the macroscopic dynamics can be performed by 
comparing with the mean-field equations J53I54I > when rogue 
dissociation is negligible (see below). As for the microscopic 
dynamics, we have to ensure that the mean-field scattering 
length 1421 . deduced from the integral expression 14 H . is com- 
puted accurately enough at each time step. To that purpose, 
we have solved the stationary version of the two-body equa- 
tions d36l > and extracted the phase shift of the continuum wave 
function in the open channel from its asymptotic behaviour. 
The scattering length can be deduced by extrapolation to zero 
collision energy. In Fig. |4j we compare the scattering length 
computed by this method, and the mean-field scattering length 
given by Eq. 1411 . One can see from the excellent agreement 
that the latter is indeed accurately computed in the stationary 
case. 



VI. APPLICATION: PHOTOASSOCIATION IN A SODIUM 
CONDENSATE 

A. Experimental conditions of the NIST experiment 

The chosen example starts from the photoassociation exper- 
iment in NIST EUl . where a condensate of sodium atoms is il- 
luminated by a cw laser of intensity / varying from 0. 14 to 1 .20 
kW/ cm 2 . The peak density is about 4 x 10 14 at/cm 3 . The fre- 
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Figure 4: Optically induced Feshbach resonance: variation of the 
sodium scattering length as a function of the detuning of the photoas- 
sociation laser in the vicinity of the v=135, 7=1 resonance |41], for a 
laser intensity of 15 kW/cm 2 , and a variable detuning A relative to the 
Di atomic resonance line. Spontaneous emission is not included in 
these calculations. Solid line: scattering length determined by phase 
shift from two-atom calculations. Dots: mean-field scattering length 
as defined by Eqs. 1411421 . calculated from the stationary states in the 
reduced pair wave approximation, see Eq. J66I67I . 



quency ( 16 913.37 cm" ) is chosen at resonance with the 7=1, 
v= 1 35 vibrational level in the 0+ (35 + 3Pi 12) potential curve of 
Na2. The detuning of the photoassociation laser relative to the 
Di resonance line is rather large (43 cm" 1 ), and corresponds to 
the binding energy of the v=135 level. The molecules formed in 
the excited electronic state have a natural width y = 27t x 18.38 
MHz due to spontaneous emission. 

Previous estimations |41] based on criterion \62\ |54] pre- 
dicted that these experimental conditions should lead to rogue 

dissociation, as (27t 2 pc^ 3 ) 1 3.2 for / =1.2 kW/cm 2 . How- 
ever, this criterion is relevant only in the coherent regime. Ac- 
cording to the analysis of section HVl the experiment is actually 
in the adiabatic regime (see Fig. |2j, as 2i2/y^ 0.04. There- 
fore, the relevant criterion is (1631 . As (27I 2 p^q) 1 « 0.03 
for I = 1.2 kW/cm 2 , one concludes that no rogue dissocia- 
tion should occur and the mean-field approximation should be 
valid, which is consistent with the experimental observation. 
However, if the intensity can be increased up to 10 kW/cm 2 , 

then (27t 2 p^) w 0.6 and rogue dissociation should occur 
(see Fig. |2j. We therefore present numerical calculations sim- 
ulating the experiment in this high-intensity regime to see the 
effects of rogue dissociation. Calculations in similar conditions 
were done by T. Gasenzer 15011 using the first-order cumulant 
equations with separable potentials. 



B. Parameters of the calculations 

For these calculations, the potentials curves for the ground 
state of Na2 X'E+ and 0+(3s+3pi/2) have been taken from 
Ref. l58ll . The ground state potential is chosen so that the scat- 
tering length is 54.9 ao, as in T. Gasenzer's work ll5(ill . Hy- 
perfine structure is not included, and the last bound level has 
a binding energy of 317.78 MHz ( ~ 10" 2 cm" 1 ). In the ex- 



cited curve, the level v=135 is bound by 49.23 cm" , the two 
neighbouring levels being v=134 at -52.95 cm" 1 and v=136 at 
-45.757 cm" 1 . The spontaneous emission term J47i for this 
excited state is set to y = 2n x 18.36 MHz to match the experi- 
ment fidll . 

The dipole moment matrix element is taken as D=2, assum- 
ing linear polarisation of the laser, neglecting r-dependence, 
and deducing the atomic lifetime from the long range coeffi- 
cient Ci= 6. 128 of Marinescu and Dalgarno [59]. The coupling 
term W in the equations d36i is therefore linked to the intensity 
lbyW = V2//(ce ). 

We model the cw laser as follows: the intensity starts from 
zero and is turned on linearly (Ib(t) = I X t/T; t < T) to be- 
come constant and equal to / after t = T = 0.5/js. This value of 
T, previously used in the calculations by T. Gasenzer 15 OIL cor- 
responds to the typical rise and fall-off times of the laser in the 
experiment of Ref. lEdll . Starting from a pure sodium conden- 
sate, with a density p(f = 0)=p g (f = 0), and therefore assuming 
p'(f = 0)=p e (f = 0) = 0, we have solved the coupled equations 
(1391401 and gU |46j as well as the mean-field equations J53I 
I54> . From these calculations, we obtain the time-variation of 
the relative number of condensate atoms 

a g (0=p g (0/p(0);a g (0) = l, (68) 
of non-condensate atoms in the open channel, 

a' g (f) = p' g (f)/p(0);a' g (0)=0 (69) 
and the relative number of atoms in the closed channel 

OCeW = Pe(0/P(0) = 2p e (f); p e (0) = 0, (70) 

which is twice the relative number of photoassociated 
molecules p e (f )■ The remaining fraction 1 — a g — a g — a e cor- 
responds to photoassociated molecules that have been deex- 
cited by spontaneous emission, yielding either cold molecules 
in the ground state or pairs of "hot atoms" that usually leave 
the trap. 

C. Results for high intensities 

We first performed calculations for laser intensities / in the 
range of the experiment (from 140 to 1200 W/cm 2 ). The 
photoassociation dynamics on resonance is reported in Fig. [5] 
showing the condensate, noncondensate fractions, as well as 
the mean-field scattering length as a function of time. Here we 
are far from the rogue dissociation limit (see dots a and b in Fig. 

As expected, the noncondensate fraction remains negligible 
and all approximations agree with the mean-field approxima- 
tion. 

Since we are in the adiabatic regime (see Fig. |2] again), 
the mean-field approximation is consistent with the usual 
two-body theory. This can be seen in the left column of 
Fig- <03' where we plotted the photoassociation line shape, as 
well as the variation of the mean-field scattering length as a 
function of the frequency of the laser. The line shape presents 
no difference from the one predicted by two-body theory, 
nor does the mean-field scattering length from the optically 
modified scattering length of the two-body theory flaM7ll. 
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Figure 5: On-resonance photoassociation dynamics for the v=135 level of Na2 A'Lj" using a cw laser of intensity (for t > 0.5/js) /= 140 W/cm 2 
(left column) and 7=1200 W/cm 2 (right column). This corresponds to the adiabatic regime where all approximations reduce to the mean-field 
aproximation and are consistent with usual two-body theory. Top panels (a) and (d): the decay of the condensate fraction 0C g (f) is caused by 
spontaneous emission and dissociation of the photoassociated molecules into pairs of atoms which escape the trap. This decay is more rapid in 
panel (d), due to the larger intensity (1200 W/cm 2 ) of the photoassociation laser. Middle panels (b) and (e): time variation of the fraction of 
non-condensate ground state atoms ofJt), which remains negligible. Bottom panels (c) and (f): time variation of the real 9t[a(f)] and imaginary 
3[a(f)] part of the mean-field scattering length. 



We then performed calculations for 7=10 kW/cm 2 . We are 
now in the rogue dissociation regime where deviations from 
the mean-field approximation are expected. The on-resonance 
dynamics is ploted in Fig. [6] We observe the following effects 
in both the FOC and RPW approximations: 

1. First, there is a significant final fractions of non- 
condensate atoms, which of course are ignored in the 
Gross-Pitaevskit picture. 

2. Second, the decay of the condensate at short times is 
slower than the one predicted by the Gross-Pitaevskit 
coupled equations. 

These two effects arise from the dynamic correlation <&g y " 
discussed in section HVl but have different interpretations. 

1. Creation of correlated pairs 

Let us first study how the non-condensate atoms are pro- 
duced. This can be done by analysing the ground state compo- 



nent of the pair wave function, illustrated in Fig.0 We see that 
a strong maximum emerges in the wave function <J> ? around 
55 ao, which we identify with the presence of weakly bound 
molecules, corresponding to population of the last vibrational 
levels in the X potential. However, these bound states 
give a very small contribution to the fraction of non-condensate 
atoms and they disappear when the laser is turned off. They are 
in fact a near-resonance feature of the adiabatic correlation <t>^' 
which is already explained by the stationary two-body theory. 

The significant fraction of non-condensate atoms is ex- 
plained by the appearance of waves in the pair wave function 
at larger distances (Fig.0. In both FOC and RPW calculation, 
these waves are created at short distances as soon as the laser is 
turned on and then propagate in the outward direction, reaching 
distances where they are no longer affected by the laser cou- 
pling -they are not affected when the laser is switched off. This 
explains why the fraction of non-condensate atoms becomes 
mostly constant after 5 ^s. Thus, these non-condensate atoms 
correspond to correlated pairs of free atoms with opposite mo- 
menta. Their typical kinetic energy is expected to be of the 
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Figure 6: Same as Fig. [5]for an intensity of 10 kW/cm 2 . Left column corresponds to a rise time of the laser of 0.5 fjs, whereas right column 
corresponds to a rise time of 16 jjs (the end of the ramp is indicated by a dotted-dashed vertical line). 



order of HT ~ 2 /jK, which is consistent with the wave length 
of the waves observed in the pair wave function. For weak 
trapping potentials, they may leave the system. Otherwise, one 
should treat the collisions between these correlated pairs and 
the condensate atoms in the presence of the laser field. How- 
ever this goes beyond both the FOC and RPW approximations. 



The appearance of the correlated pairs at high intensity has 
already been predicted in T. Gasenzer's work and it was shown 
that it leads to a saturation of the number of possible stable 
molecules formed by spontaneous emission 1 50] . However, it 
is worth noting that the formation of these pairs is a conse- 
quence of the fast rise (T=0.5 /js) of the laser intensity, and is 
most probably related to a two-body non-adiabatic effect - the 
waves also appear in a time-dependent two-body calculation. 
The right column of Fig. [5] shows the on-resonance dynamics 
when the laser intensity is rised more slowly for T = 16 fjs. 
We can see that the final fraction of non-condensate atoms is 
notably reduced in both approximations. This shows that the 
creation of correlated pairs is greatly sensitive to the way the 
laser is turned on. In contrast, the depletion of the condensate 
still shows a marked deviation from the mean-field prediction. 
We can then proceed with the discussion of this depletion, ig- 
noring the non-condensate atoms. 



2. Limitation of the decay rate 

In all models, the decay rate of the condensate K(t) = 
^4|3 [#m(0] I i s proportional to the imaginary part of the mean- 
field scattering length, the latter being one fourth of the pho- 
toassociation characteristic length introduced in Ref. |4l|]. The 
imaginary part of the mean-field scattering is plotted in the bot- 
tom panels of Fig. [6] At short times, both the FOC and RPW 
approximations lead to the same rate, which is clearly smaller 
than the one predicted by the mean-field approximation. Such 
a limitation of the rate was already pointed out in the coherent 
regime 15411 . and also observed in the adiabatic regime 1 50 ] . 

Our interpretation is the following: because of the rogue 
dissociation, excited molecules are coupled back to the con- 
densate through the dynamic correlation. This extra correla- 
tion reduces the efficiency of the coupling between the con- 
densate and the excited molecules, thereby reducing the loss 
from the condensate. The fact that this effect does not depend 
on the way the laser is turned on suggests that it is a genuine 
many-body effect, namely the influence of the dynamics of the 
medium on the pair dynamics. We also note that the dynamic 
correlation responsible for this limitation of the rate gives a 
very small contribution to the non-condensate atoms, unlike 
the part corresponding to the correlated free atoms. We suspect 
from Eqs. j49!52t that the Fourier transform of this correlation 
scales as 1/k 4 . 
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Figure 7: Modulus of the two components <|)g(r) and § e (r) of the pair 
wave function in the FOC approximation, as a function of the inter- 
atomic distance. Note that a logarithmic scale is used on both axes. 
Upper graph: before applying the laser: only ground state pairs ex- 
ist, § e (r) = 0; the wave function § g (r) corresponds to a zero-energy 
scattering state: it oscillates in the short distance region, where the 
potential U g (r) is not negligible; for distances r >1000 ao, the am- 
plitude is equal to the initial condensate density prj. Lower graph: 
after applying a cw laser of intensity /= 10 kW/cm 2 during 40 fis, a 
level of the excited potential has been populated: <j> e (r) oscillates up 
to the classical turning point of this level; the population transfer to 
this excited state is visible in the decrease of both the short-range am- 
plitude of § g (r), and its very long-range amplitude corresponding to 
the condensate density. Note that a strong maximum emerges around 
r ~ 55ay; it is the near-resonance signature of a weakly bound state in 
the ground-state channel. At large distances, § g (r) shows oscillations, 
interpreted as outward motion of "hot" correlated pairs (T ~ 4 to 40 



Note that at longer times in the reduced pair wave approxi- 
mation, the rate (and more generally the mean-field scattering 
length) approaches the standard two-body rate (standard scat- 
tering length A). On the other hand, in the first-order cumulant 
approximation the rate finally exceeds the standard rate and 
starts oscillating. However, these differences are not so signifi- 
cant from the experimental point of view because they occur at 
a point where there are very few condensate atoms remaining 
in the system. 



3. Symmetry of the line shape 

The differences between the RPW anf FOC approximations 
are more conspicuous in the line shape of the resonance - see 
right column of Fig. [8] In the RPW approximation the final 
condensate and non-condensate fractions are symmetric with 
respect to the resonance. On the other hand, in the FOC ap- 
proximation the line shape becomes very asymmetric: on the 
"blue" side of the resonance (where the bare resonant state lies 
below the dressed threshold) the condensate is more depleted 
and more non-condensate atoms are produced than on the "red" 
side. We suspect that this asymmetry is due to the strong de- 



pendence of the scattering properties on the mean-field energy 
in the FOC approximation, as we saw in Section lTl D 21 Indeed, 
at high-intensity the mean-field energy becomes large and pos- 
itive on the blue side while it is large and negative on the red 
side. On the other hand, the scattering properties in the RPW 
approximation correspond mainly to those of the usual zero 
energy scattering problem. 

Comparison of these predictions with experiment may not be 
straightforward. First, there are important experimental issues 
such as inhomogeneous broadenening to overcome at these 
high intensities 14111 . Second, collisions between condensate 
and non-condensate atoms are neglected in both approxima- 
tions, and as we already pointed out in Section UTeI the RPW 
approximation brings some correction to the FOC approxima- 
tion, but only in an incomplete way. However, at a qualitative 
level, the symmetry of the experimental line shapes at these 
high intensities could show the relevance of this mean-field 
correction occurring at short interatomic distance. 



4. Realistic potentials and cw lasers 

Finally, we should note that our calculations using the FOC 
approximation with realistic potentials agree with T. Gasen- 
zer's calculations which use a separable potential 1 50], except 
for times when the condensate fraction becomes negligible. 
This shows that no major improvement is brought by the details 
of the potential in these cw laser photoassociation calculations. 
Using the simplified equations (1491521 in this case would lead 
to similar results. In fact, since we have addressed only the 
adiabatic regime (see Fig.|2j, only two equations would be suf- 
ficient, namely ( 149 \ and d5 H or ( 149 \ and i52i , as the molecular 
condensate ^PmCan be eliminated adiabatically in these equa- 
tions. Note that these equations are very similar to those of 
say, Ref. [54], but are free of any ultraviolet divergence or 
renormalization process as the adiabatic correlation have been 
eliminated from the equations in the first place. 

In the more complex photoassociation processes involving 
pulsed lasers, the simplified equations are inadequate and we 
shall use the more general equations and numerical procedures 
described in this paper to investigate these situations. 



VII. CONCLUSION 

We have compared two many-body approximations for the 
time-dependent description of photoassociation and optically- 
induced Feshbach resonances in an atomic condensate: the 
first-order cumulant approximation and the reduced pair wave 
approximation. These approximations differ only in the way 
the influence of the mean field on a pair of condensate atoms is 
treated at short separations. Each approximation leads to a set 
of coupled equations describing the microscopic and macro- 
scopic dynamics for a two-channel problem. We demonstrated 
that these equations can be solved numerically with realistic 
molecular potentials, which should prove essential when ad- 
dressing experiments with chirped laser pulses, for which de- 
tails of the interaction potentials may influence the dynamics. 

From the general equations, we identified several regimes. 
We define the adiabatic regime when the excited molecular 
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Figure 8: Photoassociation line shape and optical Feshbach resonance: condensate, non-condensate fractions and mean-field scattering length 
after 40 /js, as a function of the laser detuning 8 from the bare molecular state (see Fig.0. At high intensity (10 kW/cm 2 ), the first-order cumulant 
equations predict highly asymmetric behaviours with respect to the resonance. 



channel is just an intermediate state during the collision of two 
condensate atoms. In contrast, we define the coherent regime 
when the excited channel gives rise to a molecular condensate 
having the features of a one-body condensate. In each of these 
two regimes, a mean-field theory is obtained when the pair 
correlation in the ground channel can be eliminated adiabati- 
cally. In the coherent regime, this leads to two coupled Gross- 
Pitaevskit equations B25II . In the adiabatic regime, this leads 
to a single Gross-Pitaevskit equation with a complex scattering 
length predicted by two-body theories. This is the usual regime 
investigated so far experimentally 1 20, 41, 55]. 

The condition for the breakdown of the mean-field approx- 
imation (or so-called rogue dissociation J3) is different in 
each reg ime. We showed that, contrary to previous estimates 
l4ll l54ll . the condition for current experiments are under the 
adiabatic regime. 

Solving the general equations numerically, we investigated 
the case of a sodium condensate where a photoassociation cw 
laser is turned on in conditions similar to the experiment of 
McKenzie et al |4l]]. At h igh intensities (^10 kW/cm 2 ), which 
could not be attained in 14 ill and where rogue dissociation is 
expected to occur, we observed the following effects: 

1 . Ground-state correlated pairs of free atoms are produced. 
This is agrees with the saturation predicted by T. Gasen- 
zer 15011 . The creation of such atom pairs can indeed 
strongly limit the final yield of ground-state molecules 
formed by spontaneous decay. However, Ref. 1 50] con- 
sidered only a rapid turn on of the laser. By introducing a 



slower switching procedure, leading to a more adiabatic 
behaviour, the present work shows that the production of 
hot atom pairs may be decreased. 

2. The photoassociation rate of the condensate at short 
times is smaller than the rate predicted in the mean-field 
approximation. This limitation happens at much lower 
intensities than the unitary limit of the two-body theory. 
This effect has already been predicted by J. Javanainen 
and M. Mackie 1 54] in the coherent regime. We interpret 
it as the appearance of a dynamic many-body correlation 
in the condensate. 

3. The photoassociation line shape becomes asymmetric in 
the first-order cumulant approximation, while it remains 
symmetric in the reduced pair wave approximation. 

We think this last point should be a good qualitative effect to 
experimentally distinguish between the two approximations. 
Future work will investigate the time dependence of the 
formation of stable molecules, via two-colour or one-colour 
experiments with pulsed lasers, using the theoretical and 
numerical tools developed in the present paper. 

Acknowledgements P.N. is very grateful to Eite Tiesinga for 
his essential remarks on this work. Discussions with Arkady 
Shanenko, Paul Julienne, Thorsten Kohler, Christiane Koch, 
Eliane Luc-Koenig, Philippe Pellegrini are gratefully acknowl- 
edged as well. This work was performed in the framework of 



17 



the European Research Training Network "Cold Molecules", 
funded by the European Commission under contract HPRN CT 
2002 00290. RN. acknowledges for an invitation in Clarendon 
Laboratory, Oxford. 



VIII. APPENDIX 

A. Definition of the non-commutative cumulants 

If A i , A 2 , • • A„ are n bosonic operators, the n-order cumulant 
{A\A 2 ..A n ) c is defined recursively as follows [34]: 

(Ai) = {AiY 
(A X A 2 ) = (A l r(A 2 r + (A l A 2 y 
{A X A 2 A 3 ) = (A l ) c (A 2 y\A 3 ) c + (AtfiA^y 
+ (A 2 ) 6 '(AiA 3 ) f + (A 3 ) c (AiA 2 ) t 
+ (AiA 2 A 3 ) c 



where (...) denotes the quantum average in the many -body 
states of the system. For an ideal gas, we can check that all cu- 
mulants containing more than two bosonic operators are zero, 
according to Wick's theorem. For a system close to the ideal 
gas, the cumulants are non-zero but tend to zero as the order is 
increased. They are therefore a sort of measure of the deviation 
from the ideal case. 



B. In-medium effective wave functions 

1. General expressions 

The in-medium pair wave function approach, as we may call 
it, was initiated by a series of papers by A. Yu. Cherny and A. 
A. Shanenko. Its aim is to have a many-body description of the 
dilute Bose gas which remains valid at short interatomic dis- 
tances, so that interaction potentials with strong repulsive cores 
can be treated directly. The authors have mainly addressed 
the problem of the ground state for a homogeneous system. 
In Ref. 13 811 . we have generalized some of their ideas to the 
inhomogeneous time-dependent case. We will recall here the 
derivation of the in-medium pair wave functions given in 1 38], 
and will show in addition how three-body wave functions can 
be constructed out of these pair wave functions. This will lead 
to the reduced pair wave approximation used in this paper. 

The starting point is to consider the reduced density matrices 
of the the many -boson system. As these matrices are hermitian, 
they can be diagonalized in a basis of orthogonal eigenvectors, 
associated with positive eigenvalues. For example, the follow- 
ing one-body, and two-body reduced density matrices can be 
diagonalized as follows: 

W + (x)<j/(y)) = £^(x)¥,-(y) (71) 

i 

l(yV)y + (z)y(x)y(y)) = Y>*(w,z)4>,(x,y) (72) 

i 

where the 'P, are the eigenvectors of the one-body density 
matrix, and 4>, are the eigenvectors of the two-body density 
matrix. We have normalized these vectors to their respective 
eigenvalues, which means thatM = /«i 3 x| l P,(x)| 2 isthe eigen- 
value associated with I*,-, and M; = / c/ 3 yc/ 3 x|<I>,(x,y)| 2 is the 
eigenvalue associated with <!>,. As these eigenvectors have the 
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form of wave functions, we call them effective one-body wave 
functions and effective pair wave functions. Nj is interpreted as 
the average number of particles of the system in the one-body 
state |\Pj), and M, is interpreted as the average number of pairs 
in the two-body state |<t>,). We can check that = N and 
Y,iMj = jN(N— 1), where N is the total number of particles in 
the system. 

In the case of a condensate, a certain one-body wave func- 
tion ^Po has a macroscopic occupation number No S> L/^o^'- 
In the 1/(1) symmetry breaking picture, this one-body wave 
function is identified with the condensate wave function, or or- 
der parameter (Vj/) . According to the Bogoliubov prescription, 
the field operator Vff can then be decomposed into its average 
value and a remaining fluctuating operator 8. Expanding 
the two-body density matrix i72l with the Bogoliubov pre- 
scription, and refactorising the expression, we showed that it 
can be written: 



L*w(w,i)*w(x,y) 



(73) 



i=0 



-i(0 t (w)§ t (z)§(x)e(y))-f:< I (w,zX(x,y) 



where: 



^oo(x,y) 

<%( x -y) 

for i 7^ and: 



— ^o(x)^ (y) + <&'oo(x,y) 
t P (x)^(y)+ t Po(y)^(x) 

V2 



(74) 



<(x,y) (75) 



<& 00 (x, y ) - — (9(x)9(x)) 



<(x,y) 



d 3 z^4(e + (z)6(x)e(y)) for/^0 



As the 4>o,'s and <f>oo are orthogonal in the limit of large 
systems, we conclude that we have found the first pair wave 
functions of the diagonalized form (I72> . 

<Poo is interpreted as the condensate-condensate wave func- 
tion, ie the wave function for a pair of particles both com- 
ing from the condensed part. It is composed of two terms: 
an asymptotic term 4-\Po(x)\Po(y) corresponding to the free 
motion of two independent condensate particles, and a scat- 
tering term <t> 00 corresponding to their correlated motion due 
to their interaction. Similarly, <t>o; for i ^ is interpreted 
as the condensate-noncondensate wave function, ie the wave 
function for two particles, one coming from the condensed 
part and the other coming from the non-condensed part. It 
is also composed of an asymptotic term <p ° ( - x ^' ^J^ < " yS> ' ¥ ' ^ 

(note the symmetrization) and a scattering term <t> 0i .. We can 
check that for a sufficiently large system the norm of <t>oo is 
jA^q, the number of condensate pairs, and the norm of <t>oi is 
NqNu the number of pairs involving a condensate particle and 
a non-condensate particle in state | V F I ). The presence of a con- 
densate therefore implies a condensation of unbound pairs, as 
jN$ » LyoAW,-. 

The last line of Eq. i73l corresponds to the remaining pair 
wave functions of the system. We have no explicit expression 
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Figure 9: Schematic reprensentations of some effective in-medium 
wave functions in a condensate: the condensate pair wave func- 
tion (a), the condensate-noncondensate pair wave function (b), 
the condensate three-body wave function (c), and the condensate- 
condensate-noncondensate three-body wave function (d) . The circles 
indicate the asymptotic one-body behaviours, and the dashed lines 
represent the correlation. 



for those, but if we assume that the non-condensate particles 
form an ideal gas (ie do not interact), we can use Wick's theo- 
rem to express the last line, and check that it is equal to: 



with: 



*y(*,y) 



£ <D?-(w,z)<D y (x,y) 

0<i<j 



V2 



for i ^ j 



which are indeed the expected pair wave functions for two non- 
interacting non-condensate atoms 1 62] . 

For convenience, we may rewrite the scattering terms ap- 
pearing in i76i and 177t with multiplicative correlation func- 
tions (po/ defined as follows: 



*oo(x,y) 



1 



V2 



^o(x)^ (y)cpoo(x,y) 



_ , . «P (x)^(y)+^o(y)^/(x) , , 
<Mx,y) = — w 9o/(x,y) 



(76) 
(77) 



We call these functions cpo/ the reduced pair wave functions. 
Assuming that the particles must be decorrelated at large dis- 
tances, we expect that the reduced pair wave functionstend to 
1 at large distances. Thus, we may write cpo; = 1 + cp 0; for 
convenience, with <p 0l - vanishing at large distances. Another 
important property is that for an interaction with a strong re- 
pulsive core, the probability of finding two particles close to 
each other should tend to zero at short distances; we therefore 
expect that the reduced pair wave functionscpo; tend to zero at 
short distances. This is precisely because they tend to zero at 
short distances that they have a regularising effect when multi- 
plied with the interaction potential. 

We can now express the quantum averages found in Eqs. (|9jl 
and dl4> in terms of pair wave functions. By expanding the 
quantum average in Eq. © with the Bogoliubov prescription 
and refactorising the terms, we find: 

(<j/t( z )xpr(z)<j/(x)) = V2~f>*(z)<I>o«(z,x) (78) 

This gives a clear interpretation of this term in Eq. l|9}: the 
condensate particles can collide either with another condensate 
particle or with a non-condensate particle. Each type of col- 
lision gives rise to a term looking like a scattering amplitude, 
where the corresponding two-body wave function is involved. 
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Similarly, by expanding the quantum average in Eq. fl!4i in with the three-body wave functions: 
terms of cumulants, and refactorising the terms, we find: 

= V3lf V*(z)£l m (z,x,y) (79) 

;=o 

I 



£2ooo(z,x,y) 



£2oo;(z,x,y) = 



-^VatyWo^ofy) (l+(p'oo(z,x)+(p' 00 (z,y)+(p' 00 (x,y)J + £2' 000 (z,x,y) 



1 viWHoWPqW U+^x)+tfoi&y)+<?Uxjj) ^ 

+ ^,(y)^o(z)^o(x) M +(p^(z,x) + <(z,y) +(p' 0i (x,y)j 
v + ^,(x)^o(z)^o(y) (l +«z,x) +(p^(z,y) +(p' 0i (x,y)) ) 

I 



^ooi( z > x .y) 



where: 



£2' 000 (z,x,y) = ^=(0(z)0(x)0(y)) e 



^oo/(z,x,y) = /fiPw 



3\Ni 



e + (w)e(z)e(x)e(y)) c - 



We note that the first term of each three-body wave functions 
is not the wave function of independent particles, contrary to 
what we found with pair wave functions. Indeed, it already 
contains some pair correlation through the terms (p'. However, 
these pair correlations are not complete. For an interaction 



with a strong repulsive core, the probability for finding three 
atoms very close to one another must tend to zero. But here, 
the first term of each three-body wave function does not van- 
ish when x w y « z. We think this is because this first term is 
only the asymptotic behaviour of the three-body wave function 
at large distances. The terms £2' must contain stronger correla- 
tions at short distances. We expect that at short distances, the 
three-body wave functions are in fact fully correlated through a 
product of reduced pair wave functions akin to a Jastrow wave 
function|46|: 



£2ooo(z,x,y) = 
£2 o/(z,x,y) = 



-^=1 , o(z) , Po(x) , Po(y)9oo(z,x)(poo(z,y)(poo(x,y) + flg D0 (z,«,y) 




(z)¥ (x) l Po(y) <Poi (z,x)q>oi (z,y)<Poo (x,y) 
¥i (y) 1 !^ (z)^o (x) <poo (z,x)(p ,(z,y)cpo,(x,y) 
% (x)¥ (z)^o (y ) (z ,x) cpoo (z,y)9oi (x,y ) 



+ a' ' 0i (z,x,y) 



This is the most natural structure which preserves the sym- 
metry of the three-body wave functions and leads to a zero 
probability at short distances. The remaining terms £2" are sup- 
posed to contain the three-body correlations which cannot be 
expressed in terms of two-body correlations. Note again that 
the norm of £2ooo is ^jA^ corresponding to the number of con- 
densate triplets, while the norm of £2oo; is ^N^Nj, correspond- 
ing to the number of triplets involving two condensate particles 
and one particle in the non-condensate state |\P,-). 



The interpretation of ( I79> appearing in Eq. dl4> is again very 
simple: the pairs of condensate atoms can collide either with 
another condensate atom, or with a non-condensate atom. Each 
type of collision is accounted for by a scattering amplitude-like 
term, involving the corresponding three-body wave function. 



2. Simplified expressions 

Neglecting the collisions with non-condensate particles as 
well as the three-body correlation £2q 00 , we can express the 
quantum averages in terms of only the condensate wave func- 
tion *P = and the condensate pair wave function <t> = 
\/24> 00 : 

(\j/ t (w)\j/(z)\j/(x)) « *F(w)<I>(z,x) (80) 
(Vj/ + (w)*j/ + (z)y(x)y(y)) » 4>*(w,z)4>(x,y) (81) 

m^m^my)) « y( w) *feg*j*^ (82) 

Note that all these expressions respect the symmetry of the 
quantum averages by exchange of coordinates. Approxima- 
tions J80l > and (18 1 1 are found within the first-order cumulant 
approach, but not Approximation i82l . The reduced pair wave 
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approximation J27I used in this paper is a simplified version of fl82i where the correlation between z and x is neglected. 
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